<!DOCTYPE html>
<html>
<head>
<title>周志华《机器学习》课后习题解答系列（四）：Ch3.3 - 编程实现对率回归</title>
<meta http-equiv="Content-Type" content="text/html; charset=utf-8" />
<style type="text/css">
/* GitHub stylesheet for MarkdownPad (http://markdownpad.com) */
/* Author: Nicolas Hery - http://nicolashery.com */
/* Version: b13fe65ca28d2e568c6ed5d7f06581183df8f2ff */
/* Source: https://github.com/nicolahery/markdownpad-github */

/* RESET
=============================================================================*/

html, body, div, span, applet, object, iframe, h1, h2, h3, h4, h5, h6, p, blockquote, pre, a, abbr, acronym, address, big, cite, code, del, dfn, em, img, ins, kbd, q, s, samp, small, strike, strong, sub, sup, tt, var, b, u, i, center, dl, dt, dd, ol, ul, li, fieldset, form, label, legend, table, caption, tbody, tfoot, thead, tr, th, td, article, aside, canvas, details, embed, figure, figcaption, footer, header, hgroup, menu, nav, output, ruby, section, summary, time, mark, audio, video {
  margin: 0;
  padding: 0;
  border: 0;
}

/* BODY
=============================================================================*/

body {
  font-family: Helvetica, arial, freesans, clean, sans-serif;
  font-size: 14px;
  line-height: 1.6;
  color: #333;
  background-color: #fff;
  padding: 20px;
  max-width: 960px;
  margin: 0 auto;
}

body>*:first-child {
  margin-top: 0 !important;
}

body>*:last-child {
  margin-bottom: 0 !important;
}

/* BLOCKS
=============================================================================*/

p, blockquote, ul, ol, dl, table, pre {
  margin: 15px 0;
}

/* HEADERS
=============================================================================*/

h1, h2, h3, h4, h5, h6 {
  margin: 20px 0 10px;
  padding: 0;
  font-weight: bold;
  -webkit-font-smoothing: antialiased;
}

h1 tt, h1 code, h2 tt, h2 code, h3 tt, h3 code, h4 tt, h4 code, h5 tt, h5 code, h6 tt, h6 code {
  font-size: inherit;
}

h1 {
  font-size: 28px;
  color: #000;
}

h2 {
  font-size: 24px;
  border-bottom: 1px solid #ccc;
  color: #000;
}

h3 {
  font-size: 18px;
}

h4 {
  font-size: 16px;
}

h5 {
  font-size: 14px;
}

h6 {
  color: #777;
  font-size: 14px;
}

body>h2:first-child, body>h1:first-child, body>h1:first-child+h2, body>h3:first-child, body>h4:first-child, body>h5:first-child, body>h6:first-child {
  margin-top: 0;
  padding-top: 0;
}

a:first-child h1, a:first-child h2, a:first-child h3, a:first-child h4, a:first-child h5, a:first-child h6 {
  margin-top: 0;
  padding-top: 0;
}

h1+p, h2+p, h3+p, h4+p, h5+p, h6+p {
  margin-top: 10px;
}

/* LINKS
=============================================================================*/

a {
  color: #4183C4;
  text-decoration: none;
}

a:hover {
  text-decoration: underline;
}

/* LISTS
=============================================================================*/

ul, ol {
  padding-left: 30px;
}

ul li > :first-child, 
ol li > :first-child, 
ul li ul:first-of-type, 
ol li ol:first-of-type, 
ul li ol:first-of-type, 
ol li ul:first-of-type {
  margin-top: 0px;
}

ul ul, ul ol, ol ol, ol ul {
  margin-bottom: 0;
}

dl {
  padding: 0;
}

dl dt {
  font-size: 14px;
  font-weight: bold;
  font-style: italic;
  padding: 0;
  margin: 15px 0 5px;
}

dl dt:first-child {
  padding: 0;
}

dl dt>:first-child {
  margin-top: 0px;
}

dl dt>:last-child {
  margin-bottom: 0px;
}

dl dd {
  margin: 0 0 15px;
  padding: 0 15px;
}

dl dd>:first-child {
  margin-top: 0px;
}

dl dd>:last-child {
  margin-bottom: 0px;
}

/* CODE
=============================================================================*/

pre, code, tt {
  font-size: 12px;
  font-family: Consolas, "Liberation Mono", Courier, monospace;
}

code, tt {
  margin: 0 0px;
  padding: 0px 0px;
  white-space: nowrap;
  border: 1px solid #eaeaea;
  background-color: #f8f8f8;
  border-radius: 3px;
}

pre>code {
  margin: 0;
  padding: 0;
  white-space: pre;
  border: none;
  background: transparent;
}

pre {
  background-color: #f8f8f8;
  border: 1px solid #ccc;
  font-size: 13px;
  line-height: 19px;
  overflow: auto;
  padding: 6px 10px;
  border-radius: 3px;
}

pre code, pre tt {
  background-color: transparent;
  border: none;
}

kbd {
    -moz-border-bottom-colors: none;
    -moz-border-left-colors: none;
    -moz-border-right-colors: none;
    -moz-border-top-colors: none;
    background-color: #DDDDDD;
    background-image: linear-gradient(#F1F1F1, #DDDDDD);
    background-repeat: repeat-x;
    border-color: #DDDDDD #CCCCCC #CCCCCC #DDDDDD;
    border-image: none;
    border-radius: 2px 2px 2px 2px;
    border-style: solid;
    border-width: 1px;
    font-family: "Helvetica Neue",Helvetica,Arial,sans-serif;
    line-height: 10px;
    padding: 1px 4px;
}

/* QUOTES
=============================================================================*/

blockquote {
  border-left: 4px solid #DDD;
  padding: 0 15px;
  color: #777;
}

blockquote>:first-child {
  margin-top: 0px;
}

blockquote>:last-child {
  margin-bottom: 0px;
}

/* HORIZONTAL RULES
=============================================================================*/

hr {
  clear: both;
  margin: 15px 0;
  height: 0px;
  overflow: hidden;
  border: none;
  background: transparent;
  border-bottom: 4px solid #ddd;
  padding: 0;
}

/* TABLES
=============================================================================*/

table th {
  font-weight: bold;
}

table th, table td {
  border: 1px solid #ccc;
  padding: 6px 13px;
}

table tr {
  border-top: 1px solid #ccc;
  background-color: #fff;
}

table tr:nth-child(2n) {
  background-color: #f8f8f8;
}

/* IMAGES
=============================================================================*/

img {
  max-width: 100%
}
</style>
</head>
<body>
<p>这里采用<strong>Python-sklearn</strong>的方式，环境搭建可参考<a href="http://blog.csdn.net/snoopy_yuan/article/details/61211639"> 数据挖掘入门：Python开发环境搭建（eclipse-pydev模式）</a>.</p>
<p>相关答案和源代码托管在我的Github上：<a href="https://github.com/PY131/Machine-Learning_ZhouZhihua">PY131/Machine-Learning_ZhouZhihua</a>.</p>
<h4>3.3 编程实现对率回归</h4>
<blockquote>
<p><img src="Ch3/3.3.png" /></p>
</blockquote>
<p>所使用的数据集如下：</p>
<blockquote>
<p><img src="Ch3/3.3.1.png" /></p>
</blockquote>
<p>本题是本书的第一个编程练习，采用了自己编程实现和调用sklearn库函数两种不同的方式，详细解答和编码过程：（<a href="">查看完整代码</a>）：</p>
<h5>1.获取数据、查看数据、预处理：</h5>
<p>观察数据可知，X包含（密度、含糖量）两个变量，y为西瓜是否好瓜分类（二分），由此生成.csv数据文件，在Python中用Numpy读取数据并采用matplotlib库可视化数据：</p>
<p>样例代码：</p>
<pre><code>'''
data importion
'''
import numpy as np
import matplotlib.pyplot as plt  

# load the CSV file as a numpy matrix
dataset = np.loadtxt('../data/watermelon_3a.csv', delimiter=&quot;,&quot;)

# separate the data from the target attributes
X = dataset[:,1:3]
y = dataset[:,3]

# draw scatter diagram to show the raw data
f1 = plt.figure(1)   
plt.title('watermelon_3a')  
plt.xlabel('density')  
plt.ylabel('ratio_sugar')  
plt.scatter(X[y == 0,0], X[y == 0,1], marker = 'o', color = 'k', s=100, label = 'bad')
plt.scatter(X[y == 1,0], X[y == 1,1], marker = 'o', color = 'g', s=100, label = 'good')
plt.legend(loc = 'upper right')  
plt.show()
</code></pre>

<blockquote>
<p><img src="Ch3/3.3.2.png" /></p>
</blockquote>
<h5>2.采用sklearn逻辑回归库函数直接拟合：</h5>
<p>虽然样本量很少，这里还是先划分训练集和测试集，采用sklearn.model_selection.train_test_split()实现，然后采用sklearn.linear_model.LogisticRegression，基于训练集直接拟合出逻辑回归模型，然后在测试集上评估模型（查看混淆矩阵和F1值）。</p>
<p>样例代码：</p>
<pre><code>''' 
using sklearn lib for logistic regression
'''
from sklearn import model_selection
from sklearn.linear_model import LogisticRegression
from sklearn import metrics


# generalization of test and train set
X_train, X_test, y_train, y_test = model_selection.train_test_split(X, y, test_size=0.5, random_state=0)

# model training
log_model = LogisticRegression() 
log_model.fit(X_train, y_train) 

# model testing
y_pred = log_model.predict(X_test)

# summarize the accuracy of fitting
print(metrics.confusion_matrix(y_test, y_pred))
print(metrics.classification_report(y_test, y_pred))
</code></pre>

<p>得出混淆矩阵和相关度量（查准率（准确率）、查全率（召回率），F1值）结果如下：</p>
<pre><code>[[4 1]
 [1 3]]
             precision    recall  f1-score   support

        0.0       0.80      0.80      0.80         5
        1.0       0.75      0.75      0.75         4

avg / total       0.78      0.78      0.78         9
</code></pre>

<p>由混淆矩阵可以看到，由于样本本身数量较少，模型拟合效果一般，总体预测精度约为0.78。为提升精度，可以采用<strong>自助法</strong>进行重抽样扩充数据集，或是采用<strong>交叉验证</strong>选择最优模型。</p>
<p>下图是采用matplotlib.contourf绘制的决策区域和边界，可以看出对率回归分类器还是成功的分出了绝大多数类：</p>
<blockquote>
<p><img src="Ch3/3.3.3.png" /></p>
</blockquote>
<h5>3.自己编程实现逻辑斯蒂回归</h5>
<p>编程实现逻辑回归的主要工作是求取参数w和b（见书p59），最常用的参数估计方法是<strong>极大似然法</strong>，由于题3.1已经证得对数似然函数（见书3.27）是凸函数，存在最优解，这里考虑采用<strong>梯度下降法</strong>来迭代寻优。</p>
<p>回顾一下Sigmoid函数，即逻辑斯蒂回归分类器的基础模型：</p>
<blockquote>
<p><img src="Ch3/3.3.4.1.png" /></p>
</blockquote>
<p>目的是基于数据集求出最优参数w和b，最常采用的是极大似然法，参数的似然函数为：</p>
<blockquote>
<p><img src="Ch3/3.3.4.2.png" /></p>
</blockquote>
<p>根据书p59，最大化上式等价于最小化下式：</p>
<blockquote>
<p><img src="Ch3/3.3.4.3.png" /></p>
</blockquote>
<p>题3.2已证上式为凸函数，一定存在最小值，但按照导数为零的解析求解方式较为困难，于是考虑采用梯度下降法来求解上式最小值时对应的参数。</p>
<p>注：梯度下降法基本知识可参考书中附录p409页，也可直接采用书中p60式3.30偏导数公式。书中关于参数迭代改变式子如下：</p>
<blockquote>
<p><img src="Ch3/3.3.4.4.png" /></p>
</blockquote>
<p>对于迭代，可每次先根据(B.16)计算出梯度▽f(β)，然后由(B.17)更新得出下一步的Δβ。</p>
<p>接下来编程实现基本的梯度下降法：</p>
<p>(1)首先编程实现对象式3.27：</p>
<pre><code>def likelihood_sub(x, y, beta):
    '''
    @param x: one sample variables
    @param y: one sample label
    @param beta: the parameter vector in 3.27
    @return: the sub_log-likelihood of 3.27  
    ''' 
    return -y * np.dot(beta, x.T) + np.math.log(1 + np.math.exp(np.dot(beta, x.T)))   

def likelihood(X, y, beta):
    '''
    @param X: the sample variables matrix
    @param y: the sample label matrix
    @param beta: the parameter vector in 3.27
    @return: the log-likelihood of 3.27  
    '''
    sum = 0
    m,n = np.shape(X)  

    for i in range(m):
        sum += likelihood_sub(X[i], y[i], beta)

    return sum   
</code></pre>

<p>(2)然后基于训练集（注意x-&gt;[x,1]），给出基于3.27似然函数的<strong>定步长梯度下降</strong>法，注意这里的<strong>偏梯度</strong>实现技巧：</p>
<pre><code>'''
@param X: X is the variable matrix 
@param y: y is the label array
@return: the best parameter estimate of 3.27
'''
def gradDscent_1(X, y):  #implementation of basic gradDscent algorithms

    h = 0.1  # step length of iteration
    max_times= 500  # give the iterative times limit
    m, n = np.shape(X)

    beta = np.zeros(n)  # parameter and initial to 0
    delta_beta = np.ones(n)*h  # delta parameter and initial to h
    llh = 0
    llh_temp = 0

    for i in range(max_times):
        beta_temp = beta

        # for partial derivative 
        for j in range(n): 
            beta[j] += delta_beta[j]
            llh_tmp = likelihood(X, y, beta)
            delta_beta[j] = -h * (llh_tmp - llh) / delta_beta[j]                
            beta[j] = beta_temp[j]

        beta += delta_beta            
        llh = likelihood(X, y, beta)

    return beta
</code></pre>

<p>通过追踪参数，查看其收敛曲线，然后来调节相关参数（步长h，迭代次数max_times）。下图是在当前参数取值下的beta曲线，可以看到其收敛良好：</p>
<blockquote>
<p><img src="Ch3/3.3.4.png" /></p>
</blockquote>
<p>(3)最后建立Sigmoid预测函数，对测试集数据进预测，得到混淆矩阵如下：</p>
<pre><code>[[ 4.  1.]
 [ 1.  3.]]
</code></pre>

<p>可以看出其总体预测精度（7/9 ≈ 0.78）与调用sklearn库得出的结果相当。</p>
<p>(4)采用<strong>随机梯度下降</strong>法来优化：上面采用的是全局定步长梯度下降法（称之为<strong>批量梯度下降</strong>），这种方法在可能会面临收敛过慢和收敛曲线波动情况的同时，每次迭代需要全局计算，计算量随数据量增大而急剧增大。所以尝试采用随机梯度下降来改善参数迭代寻优过程。</p>
<p>随机梯度下降法的核心思想是<strong>增量</strong>学习：一次只用一个新样本来更新回归系数，从而形成在线流式处理。</p>
<p>同时为了加快收敛，采用<strong>变步长</strong>的策略，h随着迭代次数逐渐减小。</p>
<p>给出变步长随机梯度下降法的代码如下：</p>
<pre><code>def gradDscent_2(X, y):  #implementation of stochastic gradDscent algorithms
    '''
    @param X: X is the variable matrix 
    @param y: y is the label array
    @return: the best parameter estimate of 3.27
    '''
    import matplotlib.pyplot as plt  

    m, n = np.shape(X)
    h = 0.5  #  step length of iterator and initial
    beta = np.zeros(n)  # parameter and initial
    delta_beta = np.ones(n) * h
    llh = 0
    llh_temp = 0

    for i in range(m):
        beta_temp = beta    # for partial derivative 

        for j in range(n): 
            h = 0.5 * 1 / (1 + i + j)  # change step length of iterator 
            beta[j] += delta_beta[j]
            llh_tmp = likelihood_sub(X[i], y[i], beta)
            delta_beta[j] = -h * (llh_tmp - llh) / delta_beta[j]   
            beta[j] = beta_temp[j]  

        beta += delta_beta    
        llh = likelihood_sub(X[i], y[i], beta)

    return beta
</code></pre>

<p>得出混淆矩阵：</p>
<pre><code>[[ 3.  2.]
 [ 0.  4.]]
</code></pre>

<p>从结果看到的是：由于这里的西瓜数据集并不大，所以随机梯度下降法采用一次遍历所得的结果不太好，参数也没有完成收敛。这里只是给出随机梯度下降法的实现样例，这种方法在大数据集下相比批量梯度法应会有明显的优势。</p>
<hr />
<p>参考链接：
由于这是本书第一个编程，索引资料较多，择其重要的一些列出如下：</p>
<ul>
<li><a href="https://kukuruku.co/post/introduction-to-machine-learning-with-python-andscikit-learn/">Introduction to Machine Learning with Python and Scikit-Learn</a></li>
<li><a href="http://scikit-learn.org/stable/index.html">scikit-learn官方主页</a></li>
<li><a href="http://matplotlib.org/index.html">matplotlib官方主页</a></li>
<li><a href="http://blog.csdn.net/lilyth_lilyth/article/details/8973972">随机梯度下降（Stochastic gradient descent）和批量梯度下降（Batch gradient descent）的公式对比、实现对比</a></li>
<li><a href="http://blog.csdn.net/zouxy09/article/details/20319673">机器学习算法与Python实践之（七）逻辑回归（Logistic Regression）</a></li>
<li><a href="http://stackoverflow.com/questions/19054923/plot-decision-boundary-matplotlib">plot decision boundary matplotlib
Ask（matplotlib决策区域和边界绘制</a></li>
<li><a href="http://blog.csdn.net/abcjennifer/article/details/19848269">Python数据可视化——散点图</a></li>
</ul>

</body>
</html>
<!-- This document was created with MarkdownPad, the Markdown editor for Windows (http://markdownpad.com) -->
